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We investigate the phase diagram of the t— J Model on a triangular lattice using a Variational 
Monte-Carlo approach. We use an extended set of Gutzwiller projected fermionic trial wave- 
functions allowing for simultaneous magnetic and superconducting order parameters. We obtain 
energies at zero doping for the spin- 1/2 Heisenberg model in very good agreement with the best 
estimates. Upon electron doping (with a hopping integral t < 0) this phase is surprisingly stable 
variationally up to n m 1.4, while the d x 2_ y 2 + id xy order parameter is rather weak and disappears 
at n ~ 1.1. For hole doping however the coplanar magnetic state is almost immediately destroyed 
and d x 2_ y 2 + id xy superconductivity survives down to n m 0.8. For lower n, between 0.2 and 0.8, 
we find saturated ferromagnetism. Moreover, there is evidence for a narrow spin density wave phase 
around n ~ 0.8. Commensurate flux phases were also considered, but these turned out not to be 
competitive at finite doping. 

PACS numbers: 74.72.-h, 71.10.Fd, 74.25.Dw 
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I. INTRODUCTION 

The discovery of high-T c superconductivity, and the 
observation 2 , that strong correlations are important in 
connection with these compounds has led to a tremen- 
dous interest in understanding strongly correlated elec- 
tron physics. In particular the two simplest models for 
strongly correlated electrons, namely the Hubbard and 
t— J models, have been the subject of intensive studies. 

For example, one question of crucial interest is the 
interplay between superconductivity and antiferromag- 
netism close to the insulating phase in the t— J model. 
The ground state of this model on the square lattice 
is known to be antiferromagnetic at half-filling and one 
of the important questions is what happens upon dop- 
ing. All the approaches to these strong coupling prob- 
lems involve approximations, and it is sometimes difficult 
to distinguish the artefact due to approximations from 
the true features of the model. However, for the case 
of the square lattice, both the variational Monte-Carlo 
method (VMC)^*^ and mean-field theories^ have found 
a d-wave superconducting phase in the the t—J model. 
A wavefunction combining antiferromagnetism and su- 
perconductivity was proposed for the Hubbard and t— J 
modelsiSiiii, allowing to reconcile the variational results 
between these two models. This wavefunction allowed for 
an excellent variational energy and order parameter and 
a range of coexistence between superconductivity and 
anti-ferromagnetism was found. Further investigations 
of this class of wavefunctions has been very fruitful for 
the square lattice. This allowed to successfully compare 
to some of the experimental features with the high-T c 
cuprate o 3 ' 23 , even if of course many questions remain re- 
garding the nature of the true ground state of the system. 

The resonating valence bond (RVB) scenario proposed 
by Anderson^ was argued to be even more relevant in 



the geometry of the triangular lattice. At half-filling, the 
lattice is a frustrated magnet: the competition between 
the exchange integrals leads to unsatisfied bonds. The 
original expectation is that quantum fluctuations might 
lead to a spin-liquid behavior. However, at half-filling, it 
appears by now that the spin-1/2 triangular lattice has 
a three sublattice coplanar magnetic orderi2*i£. Quan- 
tum fluctuations are nevertheless strong, and the sub- 
lattice magnetization is strongly reduced due to these 
fluctuations. It is therefore expected that the magnetism 
is fragile and quickly destroyed by doping and that a 
strong RVB instability is present. Indeed, RVB mean 
field theories^ 2 ^ were used for the t—J model, and 
d x 2_ y 2 +id xy pairing was found over a significant range of 
doping. The same approach and questions that arise in 
the framework of the t—J model on the square lattice are 
thus very relevant in the present frustrated lattice. The 
success of the variational approach for the square lat- 
tice suggests to investigate the same class of variational 
wavefunctions for the triangular one. 

Besides its own theoretical interest, another motiva- 
tion for understanding the physics of electrons on a tri- 
angular lattice is provided by the recent discovery of 
superconductivity at low temperature in the C0O2 lay- 
ered compounds^ (Na.5CoO2.yH2O). In these systems, 
superconductivity is observed in a range of electron dop- 
ing S between 25% and 33%2&. Na<5Co0 2 .yH 2 consists 
of two dimensional C0O2 layers separated by thick in- 
sulating layers of Na + ions and H2O molecules. It is 
a triangular net of edge sharing oxygen octahedra; Co 
ions are at the center of the octahedra forming a 2D 
triangular lattice. Takada et aliSi speculated that this 
system might be viewed as a doped spin-1/2 Mott insu- 
lator. Based on LDA calculations 2 ®, a simplified single 
band t—J picture with negative t and electron doping 
was put forward^22i22. Such systems might thus be the 
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long-sought low-temperature resonating valence bond su- 
perconductor, on a lattice which was at the basis of An- 
derson's original ideas on a possible quantum spin liquid 
stateiiii. 

We propose in this paper to study the t—J model 
within the framework of the Variational Monte-Carlo 
(VMC) method, which provides a variational upper 
bounds for the ground state energy. In contrast to mean- 
field theory, it has the advantage of exactly treating the 
no double-occupancy constraint. VMC using simple RVB 
wave- functions has been used for the triangular lattice 3 ^ 
and it was found that d x 2_ y 2 +id xy superconductivity is 
stable over a large range of doping. However, in the previ- 
ous study the fact that the t—J is magnetically ordered 
at half- filling was not taken into account. We expect 
that the frustration in the triangular lattice may lead 
to a richer phase diagram and to many different insta- 
bilities. We thus propose here to study extended wave- 
functions containing at the same time magnetism, flux 
phase and RVB instabilities, in a similar spirit as for the 
square lattic o 12 ' 13 , in order to study in detail the inter- 
play between frustrated magnetism and superconductiv- 
ity. Given the non-collinear nature of the magnetic order 
parameter compared to the case of the square lattice, 
the task is however much more complicated. We thus 
present in this work a general mean-field Hamiltonian 
which takes into account the interplay between magnetic, 
RVB and flux-phase instabilities. The resulting varia- 
tional wave- function is sampled with an extended VMC, 
which uses Pfaffian updates rather than the usual de- 
terminant updates. We show that the interplay between 
the different instabilities leads to a faithful representa- 
tion of the ground-state at half-filling, and we also find 
good variational energies upon doping. To benchmark 
our wave-function, we carry out exact diagonalizations 
on a small 12 sites cluster and compare the variational 
energies and the exact ones. Finally, a commensurate 
spin density wave is considered, and is shown to be rele- 
vant for the case of hole doping. 

The outline of the paper is as follows: in Sec. [H] 
we present the model and the numerical technique. In 
Sec. IIIII we show the variational results both for the case 
of hole and electron doping. Finally Sec. \lV\ is devoted 
to the summary and conclusions. 



II. MODEL AND METHOD 

We study the t— J model on the triangular lattice de- 
fined by the Hamiltonian: 



H t -J = -t ^ ( c ia C 3," + h - C 



J^2[Si- Sj - jUiTlA (1) 



The model describes electrons hopping with an amplitude 
t, and interacting with an antiferromagnetic exchange 



term J between nearest neighbor sites (denoted 

of a triangular lattice. denotes the spin at site i, 

Si = h c \ a®a,0 c i,P ano - ® IS the vector of Pauli matri- 
ces. Ht-j is restricted to the subspace where there are 
no doubly occupied sites. In order to simplify the connec- 
tion to the Cobaltates we set t = —1 in the following and 
present the results as a function of the electron density 
n £ [0, 2], half-filling corresponding to n = 1. n > 1 cor- 
responds to a t—J model at n = 2 — n for t = 1, by virtue 
of a particle-hole transformation. In the first part of this 
section, we emphasize on the method to construct a vari- 
ational wave-function containing both superconductivity 
and non-collinear magnetism. The wave-function allows 
to consider 3-sublattice magnetism, however, since the 
latter wavefunction is restricted to a 3 site supercell, we 
briefly comment on a second simpler variational wave- 
function type, which allows to describe commensurate 
spin order. In the second part of the section, we define 
the relevant instabilities and the corresponding order pa- 
rameters. 



A. Variational wave-function 

In order to study this model we use a variational wave- 
function built out of the ground state of the following 
mean-field like Hamiltonian: 

H M f= E (- teieljc Lcja + h.c 

+ E 

+ 5^h i .S 1 --^5^n i , <7 (2) 
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Hmf contains at the same time BCS pairing (A 
{A aja i} i j), an arbitrary external magnetic field (hj), 
and arbitrary hopping phases (6fj ) , possibly spin depen- 
dent. These variational parameters are unrestricted on 
the A, B, C sites and the corresponding bonds of a 3- 
site supercell, as shown in Fig. Q We allow both singlet 
(A[ S ~ ^) and general triplet (A^- -1 ^) pairing symme- 
tries to be present. They correspond to choosing: 
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Since Hmf is quadratic in fermion operators it can 
be solved by a Bogoliubov transformation. In the most 
general case considered here, this gives rise to a 12 x 12 
eigenvalue problem, which we solve numerically. We then 
find the ground state of Hm f 
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FIG. 1: 3-site supercell of the triangular lattice. The onsite 
magnetic variational parameters can vary independently on 
each of the site A,B and C of the supercell. The BCS pairing 
as well as the flux vary independently on each of the different 
dashed bonds. 



Here a^j^.^.-j are numerical coefficients. Note that 
\iPmf) has neither a fixed number of particles due to 
the presence of pairing, nor a fixed total S z due to the 
non-collincar magnetic order. Thus in order to use it 
for the VMC study we apply to it the following projec- 
tors: Vn which projects the wave-function on a state 
with fixed number of electrons and Vs* which projects 
the wave-function on the sector with total S z = 0. Fi- 
nally we discard all configurations with doubly occupied 
sites by applying the complete Gutzwiller projector Vg. 
The wavefunction we use as an input to our variational 
study is thus: 



IVVar) = VqVs-Vn \ll>MF) 
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Although the wavefunction JSJ looks formidable, it can be 
reduced to a form suitable for VMC calculations. Using 



we find that 



(a | t/w) = P f (Q) 

Qi.j — a (k i .kj,a i ,a j ) ~ a (kj ,k t ,<jj ,at ) 
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where -P/(Q) denotes the Pfaffian of the matrix Q. Using 
this last relation, the function (jSJ) can now be evaluated 
numerically using a Monte Carlo procedure with Pfaffian 
updates, as introduced in Ref. 8. In the particular case 
where dk,l,t,1 = ak.i,i,i = and at S z = (this happens 
if the BCS pairing is of singlet type and the magnetic 
order is collinear), the Pfaffian reduces to a simple de- 
terminant, and the methods becomes equivalent to the 
standard Variational Monte-CarloiS technique. 



The above mean field Hamiltonian and wavefunction 
contain the main physical ingredients and broken sym- 
metries we want to implement in the wavefunction. In 
order to further improve the energy and allow for out of 
plane fluctuations of the magnetic order we also add a 
nearest-neighbor spin-dependent JastrowiS term to the 
wave- function: 



(8) 



where a is an additional variational parameter. Our final 
wavefunction is thus: 
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When a < the Jastrow factor favors all configurations 
which belong to the ground state manifold of a classical 
Ising antiferromagnet on the triangular lattice. Such a 
manifold is exponentially larg and this Jastrow fac- 
tor thus provides a complementary source of spin fluctu- 
ations. 

In what follows we use the wavefunction © directly for 
the VMC, but we also examine improved wavefunctions 
with respect to (0 that can be obtained by applying one 
or more Lanczos steps^i^i^: 



\lLs) = (l + \H t -j)\^ r ) 



(10) 



with optimized^ A. Since the calculation of Lanczos step 
wave functions beyond the first step is very time consum- 
ing, most of the results we will present here were obtained 
using a single Lanczos step. 

In the following, to clearly indicate which wavefunc- 
tion we use, we will denote them in the following way: 
MF I J I nLs, where MF denotes the fields present in 
the mean-field like Hamiltonian Hmf, J is present if we 
use the Jastrow factor, n Ls denotes the presence and 
the number of Lanczos steps applied on top of the bare 
wave function. 

As usual with the VMC procedure, these general 
wave functions are now used to minimize the expecta- 
tion value of the total energy {H t -j) by changing the 
variational parameters. We used a correlated measure- 
ment techniquei2iiii2i combined with parallel processing 
to smoothen the energy landscape and use a steepest- 
descent type routine to locate the minimum of energy. 
We then define the condensation energy e c of the opti- 
mal wave function as 



^var. ^Gutzwiller i 



(11) 



where eGutzwiiicr is the energy of the Gutzwiller wave 
function, i.e. the fully projected Fermi sea at zero mag- 
netization. In some cases we had to keep a small BCS 
pairing field to avoid numerical instabilities. Let us note 
that the linear size of the Q matrix is two times larger 
than in the simpler case of determinantal update VMC. 
Therefore our largest 108 sites cluster with Pfaffian up- 
dates corresponds roughly to a 200 sites cluster using 
standard updates. 



B. Commensurate order 

Since the mean-field Hamiltonian is restricted to 
a 3 site supercell, we investigated also a second class of 
mean-field Hamiltonians based on collinear commensu- 
rate structures, which are not contained in the previous 
Hamiltonian. For this type of phase, we used a simpler 
mean-field ansatz along the lines of Ref. |3 The mean- 
field Hamiltonian written in k-space is 



H S dw = E (( £fc ~~ ^ C L C ^ + f(Q> a ) C k+Q* C ko- 



E ( A 4 T c - 



h.c. 



(12) 



where k does run over the Brillouin zone of the original 
triangular lattice, tk is the dispersion of the free electron 
Hamiltonian, and is the Fourier transform of A^-. 
Depending on f(Q, er), the ground state of the Hamilto- 
nian is a commensurate charge density wave (/(Q,cr) = 
f{Q)) or a spin density wave (f(Q,cr) = af(Q)). 

At half-filling, we considered also several commensu- 
rate flux phases with 2tt x ^ flux per plaquette, using the 

Landau gauged with p £ {2, . . . , 10} and q < p. The one 
with the lowest energy was found to be the q — l,p = 4, 
as predicted theoretically^, giving an energy close to the 
simple d x 2_ y 2 + id xy wave- function. Upon doping how- 
ever the energy of such commensurate flux phases are 
rapidly much worse than the energies of our best wave 
functions. The main reason for this poor performance 
upon doping is the rather bad kinetic energy of these 
wave functions. 



C. Characterization of the encountered instabilities 

By minimizing all the variational parameters of the 
mean- field Hamiltonians @ and 112fl on a 12 and a 48 
site lattice, we find that the relevant instabilities present 
at the mean field level consist of: 

• a 120° coplanar antiferromagnetic order (AF), rep- 
resented in Fig. |21 

• a concomitant staggered spin flux phase instability 
(SFL) with: 

0i,i+a x ,a = Sa 
&i,i+a 2 ,o- = — 
0i,i+a 3 ,a = 

These bond phase factors correspond to a spin cur- 
rent in the z direction which is staggered on ele- 
mentary triangles of the triangular lattice. This 
instability follows rather naturally, since the 120° 
antiferromagnetic state itself already displays the 
same staggered spin currents (S^ x Sj ) on the near- 
est neighbor bonds (see Fig. |2J) . The effect of this 




FIG. 2: The variational parameters hi for the coplanar 120° 
antiferromagnetic order. The spins lie in the x — y plane. The 
z-component of the vector chirality (±1) on each triangular 
plaquette forms a staggered pattern. 



instability was rather small and visible only at half- 
filling. 

• a translationally invariant superconducting phase 
with d x 2_ y 2 + id xy singlet pairing symmetry (d + ), 
as well as the <L 



id x 



We have also 



-y* -^xy \d ). 

looked extensively for triplet pairing for both elec- 
tron and hole dopings and low J/\t\ < 0.4 on a 48 
site cluster, but with no success. The minimum en- 
ergy was always found for singlet pairing symmetry. 

a ferromagnetic state with partial or full polariza- 
tion (F), 

a commensurate collinear spin density waveii 
(SDW) instability with wavevector Qat = 

(7T,-7T/V3). 



D. Order parameters 

In order to characterize the phases described by the op- 
timal wave functions after projection, we have calculated 
the following observables: 

• the sublattice magnetization of the 120° coplanar 
antiferromagnetic order: 



M AI 



TV ^ 



(-i/Var|Si|?/>v 



(VVar | 1p v 



(13) 



We have checked that the projected magnetization 
has the correct 120° symmetry. To simplify the cal- 
culations, this expectation value has been sampled 
using wavefunctions without the projector Ts* ■ We 
have checked that this gives the same result as the 
correlation function lim^oo y/Si 
with the projector Vs* ■ 
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• the z component of the vector chirality (spin twist) 
on nearest neighbor bonds: 



X = 



— V 

3N ^ 



var | (Si X Si_|_ aa 

var/ 
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(14) 



We have checked that the measured vector chirality 
has the symmetry of the staggered currents derived 
from the 120° coplanar structure (c.f. Fig. 01. 

• the amplitude of the absolute value of the collincar 
magnetization in the spin density wave wave- 
function : 



M. 
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• the amplitude of the singlet superconducting order 
parameter: 
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where 



^i,a — C i,1 C i+a a ,l ' ' ' 



t J 



(16) 
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The angular dependence of the real space correla- 
tions corresponds to those of the unprojected pair- 
ing symmetry. We have also checked that the value 
of A is independent of the choice of a and (3. 



III. RESULTS AND DISCUSSION 
A. Half-filling 

We consider in this section the Heisenberg model 
(which is the limit of the t—J model at half- filling, up 
to a constant) The comparison with the large body of 
existing results for the Heisenberg model allows us to 
benchmark the quality of our wave-function. 

Let us first briefly discuss the symmetry of the vari- 
ational parameters at half-filling. The variational mag- 
netic field minimizes the energy for the two degenerate 
120° configurations. We found that the BCS pairing 
symmetry in the presence of AF order is of d + type, 
whereas the d~ , the d x 2_ y 2 and the d xy pairings have 
close but higher energies. Finally, a staggered spin flux 
variational order improves a little bit the energy. In- 
terestingly, this latter variational order is present in the 
ground state of the classical Heisenberg model. However, 
this instability was only relevant at half-filling, and the 
energy gain when S > is not significant. The various 
energies for these wavefunctions are shown in Fig. [3] The 
AF + d + + SFL/ J/Ls wave- function is thus the best ap- 
proximation, within our variational space, of the ground 
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FIG. 3: Energy per site e = 3(Si ■ Sj) of the different 
variational wave functions for the Heisenberg model versus 
the system size N~ 3/2 , with N = 36,48,108 sites. Or- 
dered by increasing condensation energy we find: d + (open 
squares), AF + d + / J (open triangles), AF + d + + SFL/ J 
(open diamonds), AF + d + /J/Ls (full triangles) and the 
AF + d + + SFL/J/Ls (full diamonds). The stars are the 
best estimates of the ground state energy available in the 
literature^. 



TABLE I: Comparison of the average energy 3 (Si • Sj) and the 
average magnetization Maf for the Heisenberg model (t—J 
model at half-filling) in different recent works for the 36 site 
cluster and the extrapolation to the thermodynamic limit. 
The energy and the sublattice magnetization are measured 
for our best wave- function (AF + d + + SFL/J/lLs) at half- 
filling. 



(3S, ■ Sj) 



Maf 



36 sites lattice 

our best wf 
Capriotti et ali- 
exact diagS*- 

ooxoo 
our best wf 
spin-wave results- 
Capriotti et alp- 



-0.543(1) 
-0.5581 
-0.5604 



-0.532(1) 
-0.540 
-0.545 



0.38 
0.406 
0.400 



0.36 
0.25 
0.21 



state of the Heisenberg model. We compare its energy 
with other estimates of the ground state energy in the 
literature (see Fig. |3|and Table 0. The mixture of AF 
and d + instabilities is improving a lot the energy, and 
our wave-function has significantly lower energies than 
the simple d + wave-function, and has energies very close 
to the best ones available. More precisely, we find in the 
thermodynamic limit an energy per site of e = —0.52 J 
for our variational wave-function. Applying one Lanczos 
step leads to a small further improvement of the energy 
to e = — 0.53J. A summary of the energies and of the 
120° magnetization for the Heisenberg model are given 
in Tableland in Fig. El 
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Inspection of these results shows that our wave- 
function has in the thermodynamic limit an energy only 
0.013 J higher than the estimates of more sophisticated, 
but restricted to the undoped case, methods^. Indeed, 
these latter methods use pure spin variational wave- 
function as a starting point, restoring quantum fluctua- 
tions with a non-variational method. Let us point out 
however that these methods are not giving an upper 
bound on the true ground state energy, so the ground 
state energy could in principle lie between this result and 
our variational one. Since our variational wavefunction 
already gives an excellent energy it would be interesting 
to check how the (non variational) methods used to im- 
prove the energy starting with a much cruder variational 
starting point would work with our variational wavefunc- 
tion and which energy it would give. We leave this point 
for future investigation however. 

Our wavefunction shows a reduced but finite magnetic 
order that survives in the triangular Heisenberg antifer- 
romagnet (THA). The 120° magnetization of our wave- 
function is reduced by the BCS pairing down to 72% of 
the classical value (see Fig. 0) which is somewhat larger 
than the spin- wave result. Thus in addition to having 
an excellent energy, our wavefunction seems to capture 
the physics of the ground state of the Heisenberg system 
correctly. Let us note that the BCS order of the wave- 
function is destroyed by the full Gutzwiller projector at 
half-filling. So, despite the presence of a variational su- 
perconducting order parameter, the system is of course 
not superconducting at half filling. Somehow the BCS 
variational parameter helps to form singlets, which re- 
duces the amplitude of the AF order. This is very similar 
to what happens for the t—J model on the square lattice: 
the inclusion of a superconducting gap decreases the en- 
ergy and decreases also the magnetization from M w 0.9 
down to M s» 0.7, which is somewhat larger than the 
best QMC estimates (M « 0.6, see Refsi&2£). Thus the 
wave-function mixing magnetism and a RVB gap seem 
to be interesting variationally, both in the square and 
triangular lattice, to restore spin fluctuations that were 
frozen in the pure classical magnetic wave-function. For 
the triangular lattice, the present work is the first at- 
tempt to reproduce the magnetic order in the THA in 
terms of a fermionic representation, which gives results 
in good agreement with other methods. The great ad- 
vantage of this approach is of course that the fermionic 
language allows to directly consider the case of hole and 
electron doping in the AF background, which is the case 
we consider in the following sections. 



B. Electron doping: n£ [1,2] 

Very few results exist away from half-filling, so in or- 
der to have a point of comparison for our variational ap- 
proach we will compare it with exact diagonalizations on 
very small clusters. Having ascertained that our wave- 
function is indeed in good agreement with the exact re- 



sults on a small cluster, we can then use it with confidence 
to describe much larger systems and extract the physics 
of the thermodynamic limit. 

Therefore, we start by comparing on a 12 site cluster 
different wave-functions with the exact-diagonalization 
results for the case of electron doping (see Fig.0}, since 
larger lattices are not readily available. Interestingly, it 
was found that even with only 2 Lanczos steps on our 
best wave-function (AF + d + /J) the energy has almost 
converged to the exact ground state energy at half-filling. 

Note that small system size is the worst possible case 
for a VMC method since the simple variational wavefunc- 
tion is not expected to reproduce well the short distance 
correlations, as we fix the long-range magnetic correla- 
tions in our variational ansatz by imposing an on-site 
magnetic field, but we do not introduce short range cor- 
rections. Variational Monte-Carlo instead focus on the 
long distance properties, which will become dominant in 
the energy as the lattice size increases. Nevertheless, the 
short range correlations contributes significantly to the 
energy on small lattices. One can thus expect on gen- 
eral grounds the energies to become increasingly good as 
the system size increases, provided that the correct long 
range order has been implemented in the wavefunction. 
The Lanczos iterations allow to correct this local struc- 
ture of the wave-function. Here we see that by changing 
this local structure our wave-function is converging very 
fast to the ground state. This is a good indication that 
even away from half filling our wavefunction is quite ef- 
ficient in capturing the physics of the system. Actually, 
the variance of the energy per site a 2 reaches its maxi- 
mum value for doping x = ^ (a 2 = 0.006), but apply- 
ing one Lanczos step reduces drastically the variance : 
(T 2 = 0.0004. At half-filling, a variance-energy plot for 
the three functions AF + d + /J/pLs (p = 0, 1,2) allows 
to extrapolate the energy at variance, and we get an 
energy per site of e = —0.61(1), which is very close to 
the exact result e = —0.6103. 

Let us now use our wavefunction to describe large sys- 
tems away from half filling. We now focus on a 108 site 
cluster, which is the largest cluster we can treat with a 
reasonable effort. We have first measured the conden- 
sation energy per site (see Fig. [5jl for different types of 
instabilities. Very interestingly, the AF j J is even better 
than a simple d x 2_ y 2 + id xy RVB state. Moreover, the 
RVB order is only weakly increasing the condensation 
energy in presence of the antiferromagnetic background 
(AF + d + /J). This is suggesting that superconductivity 
is only weakly present in the t—J model when n > 1 
which is also confirmed by the measurement of the su- 
perconducting gap (see Fig. EJ. The superconducting 
order of our best wave-function is approximately 4 times 
smaller in amplitude and in range of stability than the 
d-wave pairing in a 10 x 10 square lattice with the same 
boundary conditions. For electron doping 5 > 0.04 we 
find that the d + BCS pairing symmetry has the same en- 
ergies as the d~ one, and also as the wavefunctions with 
d x 2_ y 2 and d xy pairings. 
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FIG. 4: Condensation energy per site versus the electron dop- 
ing for a 12 site cluster for the different variational wave- 
functions. We have done exact diagonalization (ED) for a 12 
sites cluster with same periodic boundary conditions. 
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FIG. 5: Condensation energy per site e c versus electron dop- 
ing for the 108 sites lattice. We show different wave-functions 
and also the best estimate in the literature- at half-filling 
(open diamond). 
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FIG. 6: Superconducting order parameter A for a 108 site 
triangular cluster in our best wave- function (full triangles), 
in the d + wave- function (open squares). For comparison we 
show the amplitude of the d-wave gap in a 10 x 10 square 
lattice (dashed line). 
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FIG. 7: Amplitude of the 120° magnetic order Maf measured 
in the AF + d + / J for a 108 site cluster. Inset: the amplitude 
of the staggered spin-current x hi the same wave-function. 



Very strikingly, the 120° magnetic order parameter 
is surviving up to high doping 6 = 0.4, see Fig. [7| 
Long-range magnetic order at finite doping is potentially 
caused by a limitation of the VMC method, in that it 
is not possible in our calculation to model wavefunction 
with a finite correlation length, i.e. short range 120° mag- 
netic order. In our calculation, we can either totally sup- 
press the long-range 120° and get back to the Gutzwiller 
wave-function, or use the long-range 120° magnetic order 
that is highly stabilized by the potential energy. No inter- 
mediate scenario, such as incommensurate structures, is 
yet available in our calculations, but one can only expect 
the optimization of the magnetic structure to increase 
the region of stability for magnetism. We interpret this 
finding as an indication that the hole motion is not dras- 



tically modified by the presence of the non-collincar mag- 
netic structure, so that short range magnetic correlations 
will survive up to high electron doping. For the t— J z 
model on the square lattice, it is commonly understood 
that the Ising Neel order is not surviving high doping 
because of its costs in kinetic energy: whenever a hole 
wants to move in a antiferromagnetic spin background, 
it generates a ferromagnetic cloud. Therefore, good ki- 
netic energies and Neel Ising order are not compatible. 
In our case, the 3-sublattice order imposes no such con- 
straint on the kinetic energy of the holes, because of the 
120° structure. We see this fact in our wave-function 
energies: the potential energy is improved when start- 



ing from the Gutzwiller wavefunction and adding 120° 
correlations, but the kinetic energy is unchanged. Our 
best wavefunction has a better potential energy than the 
Gutzwiller wavefunction (and also than the different CFP 
wavefunction), but it also keeps the same kinetic energy. 
Therefore, this qualitatively explains why one can stabi- 
lize the 3-sublattice magnetic phase for a large set of J 
values. Finally, we note that the staggered spin-current 
pattern is also present for doping 5 = [0, 0.3] (see Fig.0. 

Interestingly, the Jastrow variational parameter a (JHJ) 
is changing sign at 5 = 0.4 when Maf ~ : for 5 < 0.4 
(6 > 0.4) the Jastrow factor favors classical Ising (ferro-) 
anti-ferromagnetic states. The VMC results show that 
the competition between the classical Ising configuration 
on the triangular lattice and the 120° order is improving 
the energy. We argue that the classical Jastrow simu- 
lates with a good approximation quantum fluctuations 
around the 120° order. Note that the Jastrow does not 
play the same role as the BCS pairing: the BCS pairing 
forms configurations of resonating singlets, and the Jas- 
trow factor forms classical Ising configurations. It is also 
worth noting that at higher doping the Jastrow parame- 
ter is leading to a small condensation energy of O.Oli for a 
large range of doping (5 — [0.4, 0.8]). It was checked that 
this gain in energy does not decrease with the size of the 
lattice and is also present for a square lattice geometry. 
We found also that for the small clusters (12 and 48 sites) 
the system was gaining a significant amount of energy 
when having a weak ferromagnetic polarization. There- 
fore, this is suggesting that the Gutzwiller wavefunction 
is not the best approximation of the ground state of the 
t— J model in the high doping limit. Nevertheless, the 
Jastrow factor does not introduce long-range correlation 
and the variational wavefunction we introduce here is still 
a Fermi-Liquid. 

Note that the variance of the energy per site a 2 reaches 
its maximum value for doping S — 0.4 for the AF + d + / J 
with a 2 = 0.0008, and applying one Lanczos step leads 
to a 2 = 0.0004. 



C. Hole doping: n G [0, 1] 

For hole doping the scenario is strikingly different. The 
120° order is weakened in a strong d + RVB background 
and disappears at doping 5 — 0.08 (see Fig.|S]and Fig.|5J. 

When superconductivity disappears, there is a first 
order transition to a commensurate spin density wave. 
No coexistence between superconductivity and the spin 
density wave was found. Then, a ferromagnetic phases 
emerges with a strong gain of condensation energy. In- 
deed the polarized states are leading to a strong gain in 
kinetic energy. This can be understood in the simple 
picture of the Stoner model, which gives a critical onsite 
repulsion related to the density of states: U^ r = l/p(ep). 
Ferromagnetism becomes favorable if ep is sitting at a 
sharp peak of pie). In the triangular lattice the tight 
binding (TB) density of states is strongly asymmetric 
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FIG. 8: Condensation energy per site e c for different wave- 
functions in a 108 site cluster. 
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FIG. 9: Amplitude of the 120° magnetic order measured in 
our best wavefunction for the 108 site cluster (left scale, full 
triangles) and the ratio of the polarization Mf on the satu- 
rated polarization M sa t in our best wavefunction (right scale, 
full circles). We show also the absolute magnetization Msdw 
for the spin density wave wavefunction (left scale, see also 
Fig.Hnj. 



and has a sharp peak at the n = 0.5 electronic density 
lying at the Van Hove singularity. Note also that the sim- 
ple t—J model of a 3 site cluster with 2 electrons shows 
that in the t > the ground state is a singlet, whereas 
the ground state is a triplet in the t < case. This shows 
that the negative sign of t with hole doping is inducing 
ferromagnetic correlations on a very small cluster. We 
find again trace of these correlations and ferromagnetic 
tendencies in the range of electronic density n£ [0.2, 0.8] 
in our 108 site lattice. Such a ferromagnetic instability 
was also predicted in Ref. |^ by comparing the energy 
of the RVB wavefunction with an analytical calculation 
of the energy of the fully polarized state. We see that 
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FIG. 10: On-site magnetization for each site of a 108 site 
lattice for the spin density wave wavefunction. Open (filled) 
circles denotes down (up) spins. The size of each circle is 
proportional to the respective amplitude of the on-site mag- 
netization. We find that the spins forms a stripe-like pattern, 
alternating ferromagnetic bonds in the a.2 direction, and anti- 
ferromagnetic bond in the two other directions. The average 
on the lattice sites of the absolute value of the local magneti- 
zation is shown in Fig. [5] 
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FIG. 11: (color online) Cartoon picture of the phase diagram 
of the t—J model we get with t < 0. Here we sketch on 
an arbitrary scale the order parameter amplitude of the 120° 
magnetic phase, the ferromagnetic phase {F), the supercon- 
ducting d x 2_ y 2 + idxy phase {SC), the commensurate spin 
density wave (SDW). Note that for electron density n > 1.04 
and n < 0.96, the energy is degenerated, within the error bars 
due to the Monte-Carlo sampling, with the pairings d x 2_ y 2, 
d xy and d x 2_ y 2 — id xy . The pitch vector of the commensu- 
rate spin density wave is Qat = (7r, — 7r/v3), and this phase 
is depicted more in details in Fig. I10I 



tributed among the lattice sites, as expected (see Fig. HUD . 
The amplitude of the on-site magnetization is shown in 
Fig. E| as a function of doping. 



minimizing the energy by changing the variational onsite 
magnetic field leads to similar results. 

Moreover, at 6 = 0.5 doping, there is a nesting of the 
Fermi surface, with three possible Q vectors. Thus, it 
is reasonable to expect that a particle-hole instability of 
corresponding pitch vector Q is stabilized close to this 
doping. We have investigated the following instabilities : 
a commensurate charge density wave, and a spin density 
wave. Interestingly, the commensurate spin density wave 
was stabilized. Indeed, we have found that the scatter- 
ing between the k and k + Q vectors introduced in the 
Hamiltonian Hsdw (particle- hole channel) allow to gain 
kinetic energy in the range of doping S = [0.15, 0.6]. For 
sake of simplicity, we have only considered the mean-field 
Hamiltonian containing one of the three possible nesting 
vectors : Qn = (ir, — ir/\/3). Finally, the phase is sta- 
bilized, when compared to the RVB and ferromagnetic 
phases, in the window S = [0.16,0.24]. Nonetheless, no 
coexistence between superconductivity and the spin den- 
sity wave was found : the energy is minimized either 
for (A fe = 0,/(Q) ^ 0), or (A fe ^ 0, f(Q) = 0) depend- 
ing on the doping, with Afc of d x 2_ y 2 + id xy symmetry 
type in the latter case. Measuring the on-site magnetiza- 
tion value, we found that the spin density wave is form- 
ing a collinear stripe-like pattern in the spins degrees of 
freedom, whereas the charge is found to be uniformly dis- 



D. Phase diagram of the model 

Based on our wavefunction we can now give the phase 
diagram for the doped system on the triangular lattice. 
The phase diagram, summarizing the various instabili- 
ties discussed in the previous sections is show in Fig. ITTI 
This phase diagram prompts for several comments. First 
one notices immediately that the competition between 
magnetism and superconductivity in this model depends 
crucially on the sign of the hopping integral (see Fig. 1 1 1 1) . 

For both hole and electron doping, the triangular lat- 
tice has a very different phase diagram from the square 
lattice one^i. In the square lattice, the AF order disap- 
pears at 6 = 0.1 and the d-wave RVB dies at 8 = 0.4 for 
the same value of J. In the triangular lattice, a similar 
stability of superconductivity exists on the hole side, but 
the electronically doped side is resolutely dominated by 
antiferromagnetic instabilities. Our results, based on an 
improved class of wavefunctions, present marked differ- 
ences with previous approximate results for the doped 
system. On the electron side, mean-field theories would 
have suggested that the long-range magnetic order state 
undergoes a first order phase transition^ into a uniform 
d + superconducting state at 5 « 3% for values of J sim- 
ilar to those considered here. A rationalization of these 
results would be that the frustration of the lattice, which 
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was from the start the motivation of RVB as a compet- 
ing state, disfavors magnetic order. Our results, where 
the Gutzwiller projection is treated exactly within the 
residual error bars due to the statistics, are in strong dis- 
agreement with this mean-field theory. Contrarily to the 
mean- field result, magnetism is dominant and the super- 
conducting order is not favored on the electronic side. In 
addition, the t—J model on the triangular lattice was ex- 
pected to have a strong and large RVB instability, since 
the coordination number is higher than on most of the 
other lattices, and naively we would expect this to pro- 
vide an easy way to form singlets. In this work we show 
that it is not the case: for electron doping the system 
is magnetic, and for hole doping the system is supercon- 
ducting, but the range of superconductivity is not ex- 
traordinary large (S < 0.16), and smaller than on the 
square lattice. 

Previous variational approaches^ were restricted to 
pure superconducting wavefunctions d x i_ y i + id xy on a 
t — t' square lattice with t = t' . In that work it was found 
that superconductivity is stabilized up to electron doping 
S w 0.24 and hole doping S ~ 0.2 for similar albeit slightly 
different values of J/t (J/t = 0.3). In our work, for the 
case of electron doping, which corresponds to the doping 
in the cobaltite experiments, our phase diagram, using 
the larger class of wavefunctions, is clearly completely 
different from this previous result, and the stabilization 
of the superconductivity in that case was clearly an arte- 
fact of the too restricted variational subspace. As show 
in Fig . 1 1 1 1 superconductivity is strongly weakened by the 
presence of 3-sublattice magnetization and is present only 
in the range of electron doping S = [0, 0.12]. On the con- 
trary, for the case of hole doping, superconductivity had 
higher energy than ferromagnetic and spin density wave 
phases for S > 0.16. We thus confirm that the previous 
results are not an artefact of their restricted variational 
subspace, and find an acceptable agreement for the phase 
diagram. However we emphasize the presence of the spin 
density wave wavefunction that was not considered in 
the mentioned work and implies a small reduction of the 
superconductivity range. 

Our calculation thus clearly prompts for a reexamina- 
tions of the arguments on the nature of superconductivity 
in a frustrated lattice. Clearly the non-collinear nature 
of the order parameter helps making the AF order much 
more stable to electron doping than initially anticipated. 
Understanding such issues is of course a very crucial and 
challenging question. Moreover, on the triangular lat- 
tice, no significant enhanced cooperative effect between 
magnetism and superconductivity seems to be observed: 
the electron doped side has a magnetic signature, and 
the hole doped side a superconducting one, but the two 
orders seem to exclude each other as much as they can, 
contrarily to what happens for the square lattice. Even 
in the parts of the phase diagram where coexistence is 
observed, coexistence between magnetism and supercon- 
ductivity in the electron doped case shows again that 
superconductivity is decreased in the presence of strong 



long range magnetic correlations. 



IV. CONCLUSION 

In this paper we have presented a variational Monte- 
Carlo study of the t—J (J/t — 0.4 and t < 0) model on 
the triangular lattice, using extended wavefunction con- 
taining both superconductivity and non-collinear mag- 
netism, as well as flux phase instabilities. The method 
we used to construct and sample the wavefunction is 
quite general and applicable to other lattices (honey- 
comb, kagome, ladders...) as well as other symmetries 
(e.g. triplet superconductivity). It thus provides a gen- 
eral framework to tackle the competition between anti- 
ferromagnetism and superconductivity in frustrated sys- 
tems. We obtained very good variational energies at half- 
filling when comparing with other more sophisticated 
methods, specialized to the half-filled case. The fermionic 
representation of our wavefunction allows to consider hole 
and electron doping. The most stable pairing corre- 
sponds to singlet pairing. We find that d x 2_ y 2 + id xy 
superconductivity is only weakly stabilized for electron 
doping in a very small window (S = [0, 0.12]) and is much 
stronger and also appears in a wider range (S = [0, 0.16]) 
in the case of hole doping. A commensurate spin density 
wave phase is leading to a gain in kinetic energy and is 
stabilized in the small window 5 = [0.16,0.24] hole dop- 
ing. Finally, ferromagnetism emerges in a wide range 
for hole doping 6 = [0.24,0.8]. Very surprisingly, the 3- 
sublattice magnetism which is present at half-filling ex- 
tends to a very wide range of electron doping 5 = [0, 0.4] 
and is suppressed very fast in the case of hole doping 
S = [0,0.08]. The large extent of 120° order for elec- 
tron doping is responsible for the suppression of super- 
conductivity. This feature was neither observed in pre- 
vious VMC calculations, nor predicted by the mean field 
theories, and prompts for a reexamination of the question 
of the stability of magnetic order on a triangular system. 

Our results show that, for electron doping, the square 
and triangular lattices behave in a very different way. For 
the square lattice, the t—J Hamiltonian finds a domain 
of stability of superconductivity and a pairing symmetry 
that is very consistent with other methods. It is thus a 
natural candidate to investigate superconducting phases 
in systems like the cuprates. For the case of the trian- 
gular lattice, the predicted phase diagram is dominated 
by antiferromagnetic instabilities, and superconductiv- 
ity, albeit slightly present, is strongly suppressed. This 
clearly indicates that, contrarily to what was suggested 
by mean-field and previous variational calculations, the 
t—J model itself is not a good starting point to tackle the 
superconductivity of the cobaltite compounds, where su- 
perconductivity is observed in the range of electron den- 
sity n = [l + |, 1 + |] . This model must be completed 
by additional ingredients to obtain a faithful description 
of the experimental system. Two missing ingredients in 
the simple t—J model could solve this discrepancy and 
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perhaps allow to obtain a superconducting instability. 
On one hand, a strong Coulomb repulsion is expected 
in this type of compound. Such a long-range interac- 
tion is not taken into account in the t— J model. Thus a 
coulomb V term should be added to get a t— J— V model. 
On the other hand, in this paper, we have used a single- 
band model as a first step to study the Co-based oxides. 
However, it is quite possible that the multi-band effect 
plays an essential role for superconductivity^. The in- 
teraction between the three bands of the compound could 
play a non trivial role in the physics of the t—J model. 
Therefore, a study of the 3-band model could also be of 



interest. Such an analysis can be done by extending the 
methods exposed in this paper to these more complicated 
models. 
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